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Abstract 

Several methods of statistical analysis are proposed and analyzed in application for a specific task 
- extraction of the structure functions from the cross sections of deep inelastic interactions of any 
type. We formulate the method based on the orthogonal weight functions and on an optimization 
procedure of errors minimization as well as methods underlying common \ 2 minimization. We an- 
alyze effectiveness of these methods usage by comparison of the statistical parameters such as bias, 
extraction variance etc., for sample deep inelastic scattering data set. 

Introduction 

Precise extraction and analysis of the nucleon's structure functions from the deep inelastic scattering 
experiments plays great role in current understanding of the particles structure and development 
of quantum field theory. Deep inelastic scattering (DIS) processes or probing internal structure of 
nucleon by a pointlike particle at small distances and high energies provide such information about 
nucleon structure. As soon as experimental facilities to investigate polarized particle interactions 
have appeared, interest moved mainly to the study of the polarized DIS that allows us to inquire into 
the structure of nucleon's spin. Phenomenologically obtained cross sections and asymmetries contain 
such information as polarized structure functions g\ and g 2 (or sometimes called g$ or g§) or partonic 
composition of nucleon's spin. To extract these data it's necessary to count on radiative corrections 
and background effects for obtaining pure cross sections (so-called unfolding of radiative smearing, see 
e.g. [HE]) for further extraction of the structure functions. This paper is devoted to some statistical 
aspects in application to extraction of the structure functions from the Born interaction parameters. 

The cross sections for neutral and charged current deep inelastic scattering both on unpolarized 
and polarized nucleon targets can be written in the following "separable" form (see e.g. [3]): 

d 2 a/dxdy = a xy = a ^ Y k (y,x, Q 2 )W k (x,Q 2 ), 

k 

where Yk(y,x,Q 2 ) are known functions, which in the Bjorken limit depend on y only and they are 
y-polynomial. Wk(x,Q 2 ) are hadronic structure functions. For massless leptons the cross section is 
parameterized by three unpolarized structure functions 

W™ pol *™ eA (x,Q 2 ) ^ Fi(x,Q 2 ), F 2 (x,Q 2 ), F 3 (x,Q 2 ) 

and five polarized functions 

W^° la,Ilzed (x,Q 2 } g!{x,Q 2 ), g 2 {x,Q 2 ), g 3 (x,Q 2 ), g :i (x,Q 2 ) y g 3 {x,Q 2 ), g 3 (x,Q 2 ). 

One can derive desired (un)polarized structure functions operating directly with asymmetries, nev- 
ertheless, as asymmetries and corresponding cross sections are being related with each other, we pro- 
pose here methods operating only with the cross sections. Procedures and methods described below 
can be easily generalized for more specific cross section expressions, e.g. extended phenomenological 
models with larger number of the structure functions or more complicated form of the Yk functions. 
These methods can be considered as an analogue of many-dimensional interpolation within a given 
hypothesis for mean values at some points of some random variable. Constraints on applicability are 
possible though. 

The structure of the given article is the following: after initial experimental (statistical) data 
being introduced we sketch orthogonal weight function method and its optimization technique for 
the simplest case of the cross section that includes two structure functions f(x) and g(x) for method- 
ological purpose (in the case of polarized particles DIS one can get such expressions by subtracting 
corresponding cross sections with opposite spin directions). Then we offer \ 2 minimization methods 
adapted to a given task. Last section compares mentioned methods by means of numerical estimation. 
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Initial experimental data 



Consider normalized on a <7o the cross section of the simplest form 

<r*v = Y+(y)f(x) + Y-(y)g(x), (1) 

where Y± reads 

Y± = 1 ± (1 - yf. (2) 
Let's suppose that we know values of cross sections a xy (in other words - a counted number of events 
in the bin of a histogram) and its errors Aa xy at some a;,y-lattice (grid) of experimental kinematical 
points. Defined grid can be either regular or irregular, so we will distinguish the following cases: 

1. Regular (rectangular) :r,y-grid 

M = X®Y, X = xi, . . . ,x m ; Y = yi,...,y„. 

2. Regular (rectangular) :r,Q 2 -grid 

M = X®Q 2 , X = xi,...,x m ; Q 2 = Q\,...,Q 2 n . 

3. Non-regular (irregular) x,y- and x,Q 2 -grid - arbitrary set of points (bins). 

Although, it's not crucial here for the lattice to be regular (besides experimenters mostly gather 
events on irregular set of bins), but at first we implement methods for regular grids and then propose 
schemes to generalize and apply them for an irregular set of points (bins) using simple interpolation 
considerations. 



Orthogonal weight function method 

Lets begin with the rectangular x, y-grid M and take the parameters a± of some preliminary chosen 
weight function u>(y;a) in such a way that the following orthogonality condition fulfils: 

]>>(y;a ± )Y ± =0, 

v 

whereas the following requirement holds: 

J>(s/;a T )y±^0. 
v 

In another way, preceding expressions mean that ui acts as a projection operator. Then one can 
extract structure function by projection in the following form: 

f{x)± A/(x) = e *(y; a-)Y + (y) £ " {y; a - ] [<Txv ± Aga>] ' 

V ! " (3) 

g(x)±Ag(x) = = Vufa^j^iAu,,]. 

T,y^(y, a +) Y -(y)y 

For example, one may choose the simplest weight functions as a u(y;a) = 1 + ay, then a± 
parameters take the following values: 

Y. y Y±{y) 

a± = "E^W (4) 

Note that in this case a± depend only on chosen data grid. 

As presented above procedure of structure function extraction implies rough estimation of the 
standard deviation (error) values in eq. (|3])- To find correctly uncertainty in the fitted structure 
functions one should adhere to the standard procedure of the variance (noted here as D) calculation 
of random variables function, which assumes here for a given x-bin and a set of n y-bins (the same 
holds for g-function) 



-C[/( cr ( a; ,J/i), • • - ,cr(x,y n ))] 



Y2^(Vi;a-)Y + (y) 



}^u)(yi;a-) 2 D[a(x,yi)]+ 

i=l 

(5) 



+ w (f*5 a -) UJ (Vf' a-)p[(r(x, yi), a(x, %■)] y/ D[a(x, yi )]D[a{x, y,)] 
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In the case of the uncorrelated data with diagonal correlation matrix p\a[x,yi),a{x, yj)] — 5y we 
specify correct estimation for deviation of the structure function values 



&f(cr(x,yi), . . .,a(x,y n )) = 



that can be larger or smaller than 



1 



]jPw({fc;a_)K|-(j/i) 



a ^2^(Vi',a-) !i D(a(x,yi)), (6) 



^^(y l ;a_)y + (y l )^ 1 



-y^,u(y t ;a-)Aa(x,yi) 



(7) 



dependently of the u>(y; a) sign. In these formulas we assume that the lattice has no uncertainty, i.e. 
_D[a±] = as given initially. 



Optimization procedure 



To minimize errors of the structure functions that are extracted from experimental data a xy ± Aa xy 
one can apply an optimization procedure of the following type. Let's introduce auxiliary functions 

A(a)=J2<"(V>a) Y +(v)> B ( a ) = $>(!««)*'- GO. 

y y 

S(x; a) = 22 u(y; a)a xy , AS(x; a) = }^ u(y; a)Aa xy . 
y y 
Two systems of equations for f(x) and g(x) and errors Af(x) and Ag(x) are 

f A(a)f{x)+B{a)g(x) = S{x;a), 
\ A(b)f(x)+B(b)g(x) = S(x;b); 

J A(a)Af(x)+B{a)Ag(x) = AS(x;a), 
\ A(b)Af(x)+B(b)Ag(x)=AS(x;b). 

To find their solutions define the following determinants: 



A(a,b) = 



A(a) B(a) 
A{b) B(b) 



and 



8Ai(x; a, b) = 



Ai(x; a, b) = 

AS(x;a) B(a) 
AS(x;b) B(b) 



S(x;a) B{a) 
S(x;b) B(b) 



A 2 (x;a,b) = 



A(a) S(x;a) 
A{b) S(x;b) 



Then we get the solution 



/(*) = ^& A/(x) = ^i# 1 g(x 



8A2(x; a, b) = 

A 2 (x; a,b) 



A(a) AS(x;a) 
A(b) AS(x-b) 



A(a,6) 



A(a,b) 



A(a,b) 



Ag{x) 



8A2(x\ a, b) 
AM) ' 



The optimal values of the parameters a and b can be found from the condition of the errors 
minimization 

]a m\w s \Af{x)\+w g \Ag{x)\], (8) 

{a,b} 

where Wf and w g - optional weight factors. 

So, optimization procedure implies determination of the optimal parameters a,k and bk for each 
experimental point Xk in order to minimize errors in this point. As a result we have the estimation 
for the mean values of the structure functions extracted at given experimental points and rough 
estimations for corresponding errors (deviations): 



6Ai(x;a,b) . A 2 (x;a,b) 
9\ x ) 



A(a,b) 



A(a,b) 



for each k = 1, 2, . . . , n. 

It should be noted that by such solution we get only approximate values for errors (likewise 
mentioned above argument about correct deviation values), nevertheless one can easily obtain correct 
values by finding minimum solution ((SI analytically and repeating formulas (jp, (|6]). 
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A difficulty arises from the fact that the method implies experimental data on a (rectangular) 
x, y-lattice which is rarely used, furthermore as a rule experimental bins chosen for analysis and 
fitting are not uniformly distributed (e.g. see kinematics in experimental reports [HE]). First of all 
we propose to apply the same scheme to regular (rectangular) x, Q 2 -lattice. The difference between 
x,y- and x, Q 2 -data consist in principle only in redefinition of the structure functions. Let's modify 
Y expressions ([2}: 



Q 



S X 



Y±(Q*) = l±(l-2-) , Y+=Q 4 , Y_=2Q 2 



,/(*) , /» 1 y /(*) 1 Y 

a x,Q 2 ~ 2 1 2~~2 + y - 

XX S Z X Z X SX 

and use the similar test weight function io(Q 2 ;a) = 1 + aQ 2 . The same orthogonality relations take 
the form of 

J2 u{Q 2 ; o±)y± = o, ^Q 2 ' a ±) + °> 

Q 2 Q 2 
a± = 



E Q2 Q 2 Y ± 
As a result one can obtain 

J2^. Q MQ 2 ;a±) = 2^]>>(Q 2 ;a ± ) + M J_ £ w( q3 ; a± )y + - M± £ w( Q^ ; a±) y_, 
<3 2 Q 2 Q 2 Q 2 



g(x) = f(x) - x 



X) ^,Q^(0 2 ; a+) = 2^ £ c(Q 2 ; a+) - ^ ]T c(Q 2 ; a+)Y_, 
Q 2 X Q 2 SX Q 2 

J2^, Q MQ 2 ;a~) = 2^^>(Q 2 ;a_) + M £ w (g 2 ; a _)y + , 

Q 2 ^ Q 2 S ^ Q 2 

where /(x) = /(x) — <?(x). Hence 

x Eq 2 o"x,q 2 ^(Q 2 ; a+) Eq 2 w (Q 2 ; a-)Y+ + xs J2q2 cr XyQ2 uj(Q 2 ; a-) Eq 2 ^(Q 2 ; a+)Y- 
J ^ x ) — — 

2 [EQ2^(Q2;a + )E Q 2^(Q 2 ;a-)Y + +x S E Q 2 ^(Q 2 ;a-)EQ 2 ^(Q 2 ;«+)y-] 

"'z 2 [Eq 2 °x,q 2 u{Q 2 ; a-) J2q 2 u(Q 2 ; a+) - Eq 2 °x,q 2 u(Q 2 ; a+) E Q2 w(Q 2 ; a_)j 

[Eq 2 ^(Q 2 ; a+) Eq2 ^(Q 2 ; 0-)?+ + xsJ2 Q2 u{q*- _) E Q 2 ^(Q 2 ; a+)y_] 

To obtain values of the absolute error one should substitute Aa x q2 instead of o x q2. 

One of the advantages of x, Q 2 -grid consists in possibility of appropriate usage of the interpolation 
through Q 2 range, contrary to the y variable range, which combines different Q 2 — sxy values in an 
unhandy way. Special optimization procedure similar to one described above by formulas is also 
applicable over here. Used scheme will be reliable if one obtains data for rectangular x, Q 2 -region 
with narrow Q 2 range to neglect existing Q 2 -dependency of structure functions (but still with at least 
two different Q 2 values). It should be noted that for non-uniform grids summing range Eq 2 may 
depend on selected x-bin and consequently a± — » a±(x) becomes a function of the given x value. 
Nevertheless, above mentioned scheme works in the same way. 

If one has x, Q 2 -lattice of the data, which cannot be grouped in the certain x-bins, but experi- 
mental points are distributed in the vicinity of certain x-values, one can eliminate these difficulties 
by interpolation methods, e.g. as briefly mentioned below. 

Interpolation 

Here we describe simplest ways to manage with data on irregular grids. One can use the following 
Taylor formula for the cross section value near given Xi value 

o"; (x,Q)=(T {x i ,Q) + —r- -^z(x - Xi)d x a J {Xi,Q) + ... 

a (x,Q ) = weighted average = — - — = r . 
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Here the additional fraction included with the purpose of normalization - but it may be omitted. 
These formulas match the case when one has data with a small dispersion in x-values. One can 
employ some external parameterization, denoted here as a llt , to fix missing experimental x points 
by simple interpolation. But total uncertainty of the extracted f(x) and g(x) values increases by 
theoretical uncertainty of these parameterizations. 

In the case of wide Q 2 range Q 2 -dependence cannot be omitted and one can treat f(x, Q 2 ) 
approximately using similar Taylor series of / for Q 2 near given Q§ value 

fix) -> r\x, Q 2 ) = f{x, Ql) + JS^g!L (Q 2 - QDdQ.f* + ..., 

or introducing some fixing factor 8ix,Q 2 ) of a pregiven form, e.g. f(x)S(x,Q 2 ), but these ideas 
require f(x) and g{x) as well as Y± and Y± to be redesignated. 



^-minimization procedure for the best a± estimations 

Below common \ 2 methods are applied to problem of structure functions extraction. First we preserve 
usage of the o± parameters (thus implicit method, requiring introducing of the orthogonal function 
u>) and modify common procedure. Let's construct the x 2 -function of random cross-sections and 
parameters a_ , a + as follows 

71 

X 2 (cr{x,y 1 ),. . . ,a(x,y„);a-,a+) = ^ [cr(x,yi) - mi(a-,a + )}D~ 1 ij {a(x,y : j) -mj(a-,a+)], 

»,j'=i 



where Dij = p[a(x, yi),a(x, yj)]^ D[a(x, yi)]D[a(x, yj)] is the covariance matrix and rn(ci-,a+) are 
the expectation values: 



^2u{yj-,a-)aix,yj) ^o;(?/ 3 ; a+)cr( 



mi(a-,a+) = Y + ( yi )^- + Y-(yi) 3 —^ . 

^2 uiyj ; a- )Y+ fa ) ^2u>(yj;a + )Y-iyj) 

2 = 1 2 = 1 

As a± are supposed to be estimated parameters only, we may neglect dependence of these mean 
values on the random cross-section values aix,yj) regarding them as a set of initial exact numbers. 

Next step is to minimize the \ 2 function. In the case of the uncorrelated data with p[a(x, yi),a(x, yj)] 
Sij we get the following estimator equations: 



n dm 

Q~-Z D ii 1 l <T ( x 'yi) " mi(a-,a+)] = 0, 

i=l 

-^2—r 1 D^ i 1 [a(x,y t ) - m l (a-,a+)} = 0. 
i=i aa+ 

Consistent solution of the system gives x 2 -estimators a±. One can check consistency calculating the 
derivatives, which should not vanish 

"X^^^A^K^S/i) ~m l (a-,a+)] + ^ (|^) D ^ =t °> 

■A d 2 m,i ! _ drrn idrrn 

-l^dK7dT D " [^^)-m i (a_,a + )]+^— D u ^- + 0. 

i—1 i — 1 

After that one should substitute obtained a± into the following expression to get f(x) and g(x): 



^2u(yj;a-)a(x,yj) ^^(yj- a+)o(x,yj) 



fix) = 3 — , 9(x) = . 

J2u(yj;a-)Y+ ( yj ) w ; a + ) Y - (yj ) 

2=1 2=1 

Apart from these values, it's necessary to get estimations for the statistical parameters such as bias, 
deviation etc. This analysis is discussed below. 
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^-minimization procedure for the best f(x) and g(x) estimations 

Here we omit a± parameter and estimate functions f(x) and g(x) explicitly (as before we refer to x 
as to fixed bin), thus changing \ 2 arguments: 

n 

X 2 (<r(x,yi), . . .,a(x,y n );f(x),g(x)) = ^ [a(x, yi ) - m i {f{x),g{x))]D' 1 ij[a(x, yj ) - rrij(f(x), g(x))], 



where D tj = p[a(x, yi),a(x, yj )\ *J D[a(x, yi)]D[a(x, yj)] is the 
are the expectation values: 



covariancc 



matrix and m(f(x),g(x)) 



mi{f(x),g(x)) = Y + ( Vi )f(x) + Y.{ Vi )g{x). 



Next step is to minimize the x 2 function. In the case of the uncorrelated data with diagonal 
correlation matrix p[a(x,yi),a(x,yj)] — Sij we get the following estimator equations: 

n 

-Y,y+{Vi)Dt\a{x, yi ) - Y + { Vi )f(x) - y-( W )s(x)] = 0, 



X>-^)A^V(*,2/0 " Y+( yi )f(x) - Y-(vi)g(x)] = 0. 



Solving this equations we get the following consistent ^-estimates f(x) and g(a 



A = 



Er=r Y+(ViW E?=i Y + { yi )Y.{ yi )D^ 
£™=i Y + { yi )Y.{ yi )D^ Y,U Y-iv^D-: 1 



/(*) = 



9{x) = ^ 



EIU Y + ( yi )D^a(x, Vi ) ELi r+Mr-G*)^ 1 
ELi ^(W)^ 1 Er=i Y + ( yi )D^a(x, Vi ) 

EIU y+feijy-CyO^ 1 EU 

One can check consistency calculating the derivatives, which should not vanish, indeed 



(10) 



Last method resembles previous optimization procedure (even in the forms of the matrices), but 
crucial point here is the absence of a parameter (it looks like Y is a weight function itself without 
any auxiliary parameter). For the case of correlated data one can refer to the same formulas using 
simple substitutions, e.g. in the matrix form: 



YYlMDj J2 Y±(y i )DT j 1 Y ± ( yj )=YlD- 1 Y ± . 

z=l i,j=l 

Considered case corresponds to the linear \ 2 approach, so using standard statistics one can get the 
following, e.g. [4]: 



£>(/) = -^ FDpT i F = ( Y-D-'Y- ) Y+D~ 
D(g) = j-^GDG T , G 



YlL> _1 Y_ ] Y^D" 1 



Y+D _1 Y + ) Y^D' 1 - ( Y^D~ 1 Y+ ) Y+D" 1 



(11) 



The similar analytical estimation of the variance holds true for previous section, but main difference 
is that one can treat values a(x,yj) either as fixed initial numbers or random ones with statistical 
uncertainty, what leads to additional complication of the calculation (i.e. m(a_,a+) holds itself 
error) and slightly increases the total deviations Af(x) and Ag(x). Main advantage of the last \ 2 
approach (least squares method in the linear case) is that it provides consistent unbiased estimator 
for / and g with the smallest in its type estimator variance, according to the Gauss-Markov theorem 

HI- 
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Discussion and numerical results 



For the purpose of numerical analysis let's take the model parameterization for the structure functions, 
e.g. in the following simplest form: 

f(x) = C f x a f(l-xff, g(x) = C g x a °(l-x) 8 s. (12) 

Then one can construct cross section a xy and its error /S.a xy at some points xi, . . . ,x n and y\ , . . . , y m , 
referring to these quantities as expectation value and deviation of some preselected probability distri- 
bution function. Also one may include additional random bias if desired. Using such initial assump- 
tions it's easy to compare numerically above listed methods. As expected from the Gauss-Markov 
theorem, numerical analysis gives both lowest values of the \ 2 (/> Q) function and the correct minimal 
variances D(f), D(g) for the last x 2 method. At same time minimum value x 2 ( a +, a -) f° r X 2 ~ a proce- 
dure equals to the minimum x 2 (/i <?) value, although presence of the random data in the mean values 
m(a_,a+) increases the errors D[f(a)], D[g(a)]. Both orthogonal method and optimized orthogonal 
method have the same larger \ 2 values for different optimal sets of a and b and for nonoptimized 
a± parameters, though they give different D(f), D(g) values. These estimations include bias to be 
calculated analytically, contrary to the linear x 2 ~/i 9 procedure without it. This causes additional 
peculiarity - one can get in this case either larger D(f) and smaller D(g) or opposite, compared to 
the x 2 methods. Detailed analysis can be carried out analytically using standard statistical methods 
after weight function definition. It should be noted that the orthogonal weight function method gives 
expectation values f{a±) and g(a±) (defined by Q) equal to optimized expectation values /(a± op t) 
and g(a±o P t) ©, although they differs from mean x 2 f an d g values (|10p. The orthogonal weight 
functions approaches with formulas <[3j) and Q can be used as a approximate estimations. The last 
mentioned least squares procedure gives reasonable unbiased estimation, and can be used for Born 
cross section analysis. It should be noted that discussed methods can be easily generalized for more 
common case and for various specific purposes. 
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